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Using a Ginzburg-Landau model, we study the phase transition behavior of compressible Ising 
systems at constant volume by varying the temperature T and the applied magnetic field h. We 
show that two phases can coexist macroscopically in equilibrium within a closed region in the T- 
h plane. It occurence is favored near tricriticality. We find a field-induced critical point, where 
the correlation length diverges, the difference of the coexisting two phases and the surface tension 
vanish, but the isothermal magnetic susceptibility does not diverge in the mean field theory. We 
also investigate phase ordering numerically. 
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I. INTRODUCTION 

Solids are under the influence of elastic constraints and 
their phase transitions are often decisively influenced by 
couplings of the order parameter and the elastic field^. 
Such elastic effects strongly depend on the nature of the 
coupling and their understanding is crucial in technol- 
ogy. In the present work, we will focus on the phase 
transition behavior of compressible ferromagnets or anti- 
ferromagnets, which has long been studied theoretically 
in the physics community^'^'^'^'^'^'^'^ i i o 1 1 1 1 1 ^ . In real ma- 
terials, the short-range spin interactions depend on the 
distances among the spins, so the spin fluctuations are 
coupled to the elastic dilation strain. In the literature 
on this problem, the main issue has been the effect of 
the elastic coupling on the critical behavior of the spin 
system. A remarkable but subtle result of the renor- 
malization group calculations^ is that the cubic elastic 
anisotropy becomes increasingly important on approach- 
ing the critical point (which is determined in the absence 
of the anisotropy). This renormalization effect should 
trigger a first order phase transition sufficiently close to 
the critical point. Simulations have been prerformed on 
compressible Ising systems and a number of numerical 
results still remain not well understoo d^°i^^'^^ . These 
theories and simulations show that the phase transition 
depends on whether the pressure or the volume is fixed. 

In this paper, we will present a mean field theory of 
compressible Ising systems at constant volume using a 
Ginzburg-Landau free energy. Our main objectives are 
to demonstrate the presence of unique two phase coex- 
istence near the tricritical point and to examine phase 
ordering after changing the temperature. Though our 
theory is a rough approximation, it will provide overall 
phase behavior for general values of the parameters. 

The organization of this paper is as follows. In Sec. 
II, we will present a model, in which the order param- 
eter and the elastic field are coupled, and eliminate the 
elastic degrees of freedom assuming the mechanical equi- 
librium condition. In Sec. Ill, we will examine the phase 
behavior in the plane of the temperature T and the order- 
ing field h. Detailed calculations will also be given on the 
susceptibility, the correlation length, and the surface ten- 



sion. The presence of a unique field-induced critical point 
will also be reported. In Sec. IV, we will numerically in- 
tegrate the time-dependent Ginzburg-Landau equation 
in two dimensions (2D). In the appendix, we will derive 
the free energy at constant pressure (or applied stress), 
where two-phase coexistence can be realized only on lines 
in the T-h plane. 



II. THEORETICAL BACKGROUND 

A. Ginzburg-Landau free energy 

We assume that a single-component order parameter 
-0 is coupled to the elastic displacement u. We set up the 
Ginzburg-Landau free energy functional F = F{^j,u\ in 
the formic, 



F 



dr 



/o + f |VVP 



/e, 



(2.1) 



where the space integral is within the system with volume 
V. The first part /o = fo{ip) depends on ip as 

fo = l^' + l^' + l^^'-hij. (2.2) 

The coefficient r depends on the temperature T as 



AoiT-To), 



(2.3) 



where Aq is a positive constant and To is the critical tem- 
perature in the absence of the elastic coupling. The other 
coefficients are treated to be independent of T. We fix 
the other field variables such as the hydrostatic pressure. 
The coefficients v and C are positive, while u can be ei- 
ther positive or negative. The h represents a magnetic or 
electric field conjugate to ip. For antiferromagnetic mate- 
rials, no uniform field conjugate to the antiferromagnetic 
order can be realized, so h — 0. We may assume h > 
without loss of generality. If /i = 0, F is invariant with 
respect to 4' The a represents the strength of the 

coupling between Tp'^ and the dilation strain, 



ei 



(2.4) 



2 



This coupling arises when the interaction among the fluc- 
tuations of tp depends on the local lattice expansion or 
contraction. 

In cubic crystals, the elastic energy density is of the 
form, 



E 



Cl2 ^ -J. 



(2.5) 



where Cn, C12, and C44 are the usual elastic moduli as- 
sumed to be constant, and e^j — iy iUj + VjUi)/2 is the 
symmetrized strain tensor. The dependence of the elas- 
tic moduli on if)'^ can be important at low temperatures, 
however. Hereafter Vi — d/dxi. The elastic stress tensor 
Gij is expressed as 

<yii = (Cii - Ci2)eii -I- Ci2ei -I- ai/)^, 

CTij = 2Ciieij {i^j)- (2-6) 

Nonvanishing ip'^ gives rise to a change in the diago- 
nal stress components. We then obtain ^ jO'ij — 
~5F /5ui, where ip is fixed in the functional derivative 
of F with respect to Ui. Note that a constant hydro- 
static pressure po can be present in the reference state, 
where the total stress tensor is pofc — <Jij ■ 



B. Elimination of elastic field at fixed volume 

The elastic field u is determined by "0 under the me- 
chanical equilibrium condition. 



0. 



(2.7) 



Furthermore, in this paper, we impose the periodic 
boundary condition on 5u in the region < x,y, z < 
V^l'^. This can be justified when the solid boundary 
is mechanically clamped. See Appendix A for the case 
of fixed applied pressure. The space averages of the 
strains then vanish; for example, (ei) = 0. Hereafter 
(•••)= / dr{- ■ ■)/¥. The following procedure of elimi- 
nating the elastic field has been derived by many authors 
in the literature in physics and engineering^'Sdi^iiiiii. 

It is convenient to use the Fourier transformation, 
u.j{r) — u^i^ exp(zfc • r), where k is the wave vector. 
Then the Fourier component of ei is expressed as 

eifc = -a<^fc/[^i2 + Cii + C44C(fc)], (2.8) 
where Lp^ is the Fourier component of the variable, 

(^(r)=V'-(V''>. (2.9) 

The space average of tp is made to vanish. The (^{k) is 
a function of the direction of the wave vector k ~ k^^k 
and is defined by 



c(fc)-^-^fc|/(i+ai), 



(2.10) 



where is the degree of cubic anisotropy, 

= (Cll — Ci2)/C44 — 2. 



(2.11) 



We have ({k) = 1 in the isotropic elasticity = 0. After 
some calculations, we may eliminate it in F to obtain the 
free energy F = F{ip} of ip only in the fornii. 



F ^ dr 



k 



The second term on the right hand side arises from the 
elastic coupling and is negative, where 

w{k) = 0?l\Cx2 + C44 + C44C(fc)]. (2.13) 
The functional derivative of F is performed to give 



(2.14) 



where /q — dfo/dtp and the Fourier transformation of ei 
is in Eq.(2.8). In equilibrium we require SF/dip — 0. 

We further simplify our free energy. In the isotropic 
elasticity, w{k) is a constant independent of k and ei — 
—oiLp/Cii. Then F is rewritten as 



F 



dr 



/o + ^|V^P-^(^^-(^^))^ 



(2.15) 



where ip is explicitly written in terms of ip and /3 is a 
positive constant defined by 



p = 2a VC; 



11. 



(2.16) 



The presence of the space average {ip'^) is a unique aspect 
arising from elasticity. 

In cubic solids with ^a < 0, w{k) is maximized along 
one of the principal crystal axes (say, along the [100] 
direction in 3D)^. If > 0, it is maximized for kj = 1/d 
for all j (say, along [111] in 3D). Let wm be the maximum 
of ^(fc) attained along these soft directions; then, 



WM 



iCa < 0), 



= a^lK + (2 - 2/d)C44] i^a > 0). (2.17) 



where K = Cu/d + Ci2(l — l/d) is the bulk modulus. 
In 2D, w{9) — w{k) is a periodic function of the angle 
9 defined by kx/k = cos9 and ky/k = sin0 with period 
7r/2, as displayed in Fig. 1. In phase ordering processes, 
the interface normals tend to be parallel to these soft 
directions, resulting in cuboidal domains^ i^^i^^i-'^^ . If the 
spatial inhomogeneity is mostly along these soft direc- 
tions except for the edge regions of the domains, the free 
energy is approximately given by Eq.(2.15) with 



(3 = 2wM- 



(2.18) 
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FIG. 1: w(6')/w(0) in 2D as a function of S/tt for ((Cn - 
Ci2)/27^,C44/if) = (0.5,1) (a), (0.3,1) (b), (1,0.5) (c), and 
(1,0.3) (d). The maximum of w(ff) is ui(0) for < and 
u)(7r/2) for ia > 0. 



C. One phase states 

We start with the free energy Eq.(2.15). If the system 
consists of a single ordered phase in equihbrium, we have 
■02 = ^02^J gj^^ homogeneous tp is determined by 

fo = {t + u^'^ + v^^)ip - h = 0, (2.19) 

where the elastic coupling disappears. The inverse sus- 
ceptibility = {dh/d^)r is given by 

=92^0/5^2 ^r + 3uV'^ + 5w7/'^ (2.20) 

We may consider the structure factor Sk of the ther- 
mal fluctuations of the Fourier component V'fc in the bulk 
region. To calculate it, we superimpose plane wave fluc- 
tuations of ip on the homogeneous average. The increase 
of the free energy in the second order yields Sk in the 
Ornstein-Zcrnike form 

Sk = l/C{k^ + K.^), (2.21) 

where k is the inverse correlation determined by 

Ck^ = d^h/dip^ -iPip^ + (3{i!^) 

= T + {Zu - 213)^'^ + ?>viP'^ . (2.22) 

In the second line, we have set {ijp') = ijp' because of the 
existence of a single phase only. Note that Cn^ in the 
second line of Eq.(2.22) is smaller than in Eq.(2.20) 
by 2/3-02. In cubic solids, k represents the inverse cor- 
relation length for the fluctuations varying in the softest 
directions. Let r take a small negative value at ft. = in 
the case m > 0; then, '02 ^ \t\/u from Eq.(2.19), leading 
to Ck^ ^ 2(1 - /3/u)|t| from Eq.(2.22). The positivity 
of is attained only for /3 < u. Obviously, the dis- 
ordered phase with -0 = is unstab le for r < 0. The 
ordered phase with 02 = —u/2v + v? /Av"^ — r (which 



is the solution of Eq.(2.19) at ft = 0) becomes unstable 
for T > Tin. In particular, as ft — ^ 0, we find 

lim Tin = -(/3^ - u^)/4w. (2.23) 



III. TWO PHASE COEXISTENCE 
A. Two phase states 
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FIG. 2: Normalized order parameters ipi/M and ipo/M and 
normalized width of the temperature window Tn,/ PM^ versus 
normalized field h/hc in two phase coexistence, where M and 
he are defined by (3.9) and (3.12), respectively. 

We show that two phases can coexist in a temperature 
window Tc ^ < T < Tc ii the parameter, 

u = u — (3, (3.1) 

is negative"*-^ and ft is smaller than a critical field ftc, 
where Tc, t^,, and ftc will be determined below. We of 
course have m < if m < or if the system undergoes a 
first order phase transition even without the elastic cou- 
pling. For < h < he the two phases are characterized 
by "0 — "00 E^nd "01 with 0i > "00 ^ 0. As ft we 
have ■00 ~* 0, while as ft —> ftc we have 0i — "00 — *■ 0. 
We will show that the space average (^02) in the free en- 
ergy (2.15) gives rise to the two phase coexistence. If it 
were neglected, we would have the usual tricritial point at 
T = u = (see the last paragraph of this subsection)ii^. 

If the volume fraction of the phase with = -01 is 
written as 0, we have 

(^2) ^0^2,^(1 _ 0)0,2^ (3.2) 

The average free energy density (/) = F/V is given by 

(/) = 0/0(^1) + (l-</')/o(^o) 

+ \m-i^o?i<l^'-^)- (3.3) 
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FIG. 3: Phase diagrams in the r-h plane with -u > for 13 /u = 
1.5, 1.4, and 1.2 (upper plate) and with m < for P/\u\ = 
1.5, 1, and 0.5 (lower plate). The t and h are scaled by ro = 
iP/v and /lo = v{\u\/vf^'^ = |6(/3/u - l)/5|^''^ftc/12, so h/ho 
is large around h ^ he for negative -u. The system is in two 
phase coexistence inside each coexistence curve (solid line), 
while it is in a one phase state outside it. Instability curve 
(dotted line) merges each coexistence curve at the critical 
point, inside which one phase states are linearly unstable. @ 



Here the interface free energy is neglected. The mini- 
mization conditions of (/) with respect to tpi and i/'o are 
given by 

mi) - - ml - V'o)V'i = 0, (3.4) 

/^(^o) + P<jy{^l - V'o)^o = 0, (3.5) 

which are equivalent to SF/6ip = at ■;/' = -01 and 4'o- 
We also minimize (/) with respect to (p to obtain 



/o(V'i 



/o(0o) + ^{^Pl 



V'^)^(20-1) = O, (3.6) 



which means that the two phases have the same free en- 
ergy density. Note that the quadratic term (c>c 0^) in (/) 
in Eq.(3.3) is positive for -0i — tAq > 0. Thus, for small 
/o(V'i) ~ /o('0o)j a minimum of (/) can be attained as 
a function of 4> in the range [0.1]. These equations may 
be solved for the simple free energy density (2.2). By 
eliminating (p we derive the equations for ipi and tJjq as 

h/v = ^PiMy^o + ^Pif/S, (3.7) 



^u/v = + + -{^Po + (3.8) 

where u is defined by Eq.(3.1). The negativity of u is 
required by Eq.(3.8). Thus ipi and ipo are independent 
of r. As /i — > 0, we have ipQ ^ and = M, where 



M = (3|u|/4i;)i/2. 
It is convenient to express V'l and as 



(3.9) 



where q satisfies 



(3.11) 



Then q/M is a dimensionless function of h/vM^, tending 
to unity as /i ^ 0. The difference ipi — ipo — {q^ — 
12h/vq'^y/^ decreases with increasing h. A field-induced 
critical! ty is attained for h — he and t — Tc, where 



he = {8/5f/^vM^/12, 
Te = 4uM^/5 - 2/?MV5. 

The critical value of the order parameter is 

= (2/5)1/2m = {3\u\/l0vy/\ 

For small positive he — h we obtain 

V;i-Vo = ^M(i-VM'^'- 



(3.12) 
(3.13) 



(3.14) 



(3.15) 



For h > he we have a unique one phase state where ip is 
determined by Eq.(2.24). In Fig. 2, we show tpi/M and 
Tpo/M versus h/he- 

Next the volume fraction of the more ordered phase (f) 
is calculated. From Eq.(3.5) it depends on r as 



4> = {Tex - t)/t^. 

This relation holds for [3 > u and h < he with 



Tex 



/3(V? - i^l). 



(3.16) 

(3.17) 
(3.18) 



In Fig. 2, the normalized window width t^/^M'^ is also 
displayed as a function of h/he- Since (p is in the range 
< (f> < 1, the two-phase coexistence is realized in the 
window region, 



T„ < T < Te: 



(3.19) 



For T below Tex the more ordered phase starts to appear, 
and is the width of the temperature window. As /i — *■ 
0, Tex and tend to the following values, 



lim Tex 

lim Tu 

h^Q 



= 3uVl6w, (3.20) 
^M^ = 3P{f3-u)/Av. (3.21) 
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On the other hand, as h he, the upper and lower 
bounds in Eq.(3.19) meet at r = and behave as t^x = 
Tc + Pipciipi - ipo) and Tex ~ = Tc - I3ipc{ipi - i>o), 
where ipi — ipo depends on he — h as m Eq. (3.15). In 
Fig. 3, we show the phase diagrams in the r-h plane for 
M > and for u < 0, separately, where the coexisting 
curves, r = Tcx and r = Tcx — t^, and the instability 
curves are displayed. The latter are obtained by setting 
Ck^ = in Eq.(2.22) using tp determined by Eq.(2.19) 
(see the discussions above Eq.(2.23)). These curves meet 
at the corresponding critical point h — h^ and t — 
given by Eqs.(3.12) and (3.13). 

The usual theory of tricriticalitjiiiii^ starts with the 
free energy density, 

for systems with short-range interactions. For this model 
a first order phase transition line^S, appears in the T-h 
plane for u < 0. (i) The line starts from the r axis {h — 0) 
at the transition point given by r = 3u^/16w where 
■0^ = 3|u|/4w in the emerging ordered phase. These val- 
ues coincide with those in Eqs.(3.20) and (3.9) in our 
elastic model, (ii) The line ends at a field-induced crit- 
ical point, where ^'^ = 3|M|/10t;, h = %v{?,\u\/\Qvf/'^/i, 
and T = %v? /2Qv. The critical values of ip and h coincide 
with those in Eqs.(3.14) and (3.12). However, the critical 
value of T is higher than that in Eq.(3.13) by 2/3Af^/5. 

B. Magnetization, susceptibility, and specific heat 





FIG. 5: Normalized susceptibility x/Xc as a function of h/hc 
and t/tu, for /3/w = 1.2, where Xc = M/hc- It is calculated 
from Eqs.(2.20) and (3.24). It increases discontinuously at the 
phase boundary from the one phase region to the two phase 
region. 
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FIG. 6: Normalized average order parameter {-(/)) /M, nor- 
malized susceptibility x/Xc, and normalized specific heat 
Cv /Cvo versus t/t^ for = 1.2 in the limit ft — > 0, where 
Cvo = TAl/2u. 



h/h. 



FIG. 4: Normalized average magnetization {tp)/M as a 
function of h/hc and r/rm for P/u = 1.2 calculated from 
Eqs.(2.19) and (3.23), where M, he, and are defined by 
Eqs.(3.9), (3.12), and (3.18), respectively 

In the two phase states in the temperature window, 
the average order parameter is given by2£ 



(3.23) 



which is continuously connected to the solution of 
Eq.(2.22) in the one phase sates outside the window 



region. See Fig. 4 for as a function of r and h 
at P/u = 1.2. The effective isothermal susceptibility 
X = (di^jj) /dh)r is calculated from 



(3-) 



where the derivatives are performed at fixed t. See Fig. 
5 for X as a function of t and h aX (3/u— 1.2. We can see 
that X is discontinuous at the boundary of the window 
region. There is no critical divergence in x at the field- 
induced criticality attained. In particular, as /i — > 0, it 
behaves as 



X^{1- 30/4 + 2wMV3/3)/(wMV3), 



(3.25) 
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where vM'^ /?> is the value of Tcx as /i ^ 0. For r > Tcx 
we have x— 1/t at ft. = 0. Figure 6 displays the behavior 
of X on the axis in the limit ft, — ^ 0. 

Next we consider the specific heat at constant volume 
Cv — — T5^(/)/9T^ (per unit volume) arising from the 
spin degrees of freedom, where ft is fixed. In the two 
phase coexistence with ft < ftc, we use Eqs.(3.3) and 
(3.16) to obtain 

Cv = TAll2(3, (3.26) 

which is independent of ft even for ft > 0. In the one 
phase region, we have Cy = TA^iP'^/{t + Zutp'^ + bvip^), 
where '0 is determined by Eq.(2.19). In particular, at 
/i = 0, Cy = for T > and Cy = TAl/2^/u^^^^ 
for T < Tcx — Tyj. In Fig. 6, we show Cy versus r at 
ft = 0. 



C. Correlation length and surface tension 



correlation lengths go to zero as 

Ko = 'ti = (4/5)«;oo(l - h/hcf'^, (3.30) 

from Eq.(3.15). If the scattering amplitude is propor- 
tional to Sk in Eq.(2.21), it grows near the critical point 
at long wavelengths. In Fig. 7, we plot ko/kqo and 
ki/kqo versus ft/ftc. It is worth noting that the inverse 
correlation length k in the one phase region also goes to 
zero at the criticality. In its vicinity, the relations (2.19) 
and (2.22) in the one phase case give 

CK^^ih-hc)/A, (3.31) 

where the term linear in r — Tc vanishes. 

We also calculate the surface tension 7 in the two 
phase coexistence. We suppose a one-dimensional inter- 
face profile ip = ipix) changing along the x direction. It 
changes from -00 at a; = —00 and to V'l at a; = 00. From 
6F/Sip = 0, we obtain 
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FIG. 7: Inverse correlation lengths kq and ki versus h/hc 
in the coexisting two phases. They are divided by kqo in 
Eq.(3.29). Normalized surface tension 7/70 is also shown, 
where 70 is in Eq.(3.36). 

Starting with the first line of Eq.(2.27), we may cal- 
culate the inverse correlation lengths, kq and ki, in the 
coexisting two phases with ip = ipQ and ipi, respectively. 
With the aid of Eqs.(3.6)-(3.8) some calculations yield 

«o = ^(V'i-^o)'(V^i+V'o)(V'i+40o), (3.27) 

1^1 = ^(V'i-V'o)'(V^i+^o)(4^i+0o). (3.28) 

As ft — > 0, we have kq — s- kqo and ni 2koo, where 

KoQ = {v/iCf'^M^ (3.29) 

is the inverse correlation length in the disordered phase 
at T = ?;M^/3 and ft = 0. As h ^ he, the inverse 



C0 = m) - /3(V' - (0'))0- (3.32) 

We integrate the above equation as 2uj — C {dip / dx)'^ , 
where uj{ip) is the grand potential, 

^ = foW - f (0' - (0'))' - Co. (3.33) 

From Eq.(3.6) the constant Co in the right hand side can 
be chosen such that to vanishes at a; = ±00 or for both 
ip = Tpo and Tpi- Some calculations yielcl^ 

to=^iip-ipo)\ip-ipif[iip + ipo + i^if + Mil (3.34) 

which turns out to be independent of r. The surface 
tension 7 is a function of ft only. It is of the form, 

dx[uj + C{di)/dxfl2] 

-00 

= / dVV2Cw(V'). (3.35) 

In the limit ft — > it becomes 

70 = lim 7 = (uC/24)i/2^/'*. (3.36) 

On the other hand, as ft — s- ftc, lo in (3.35) behaves as 
Lj = \u\{tp - Tpof{ip - ipif/2 so that 

7/70 = (32/375)(l - ft/fte)'/^ (3.37) 

which rapidly decreases near the criticality. See Fig. 7, 
where 7/70 is plotted. 
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FIG. 8; Time evolution of tp after changing r from to — 0.3ro 
at t = (upper panel) and steady profile of ip/Mo obtained 
a,t t — 10* (lower panel) for j3/u = 1.5 and h — Q in isotropic 
elasticity. Here Mq is in Eq.(4.2) and space and time are 
measured in units of I and to in Eq.(4.3). In the initial stage 
three regions with tp = ipi (black), tp = —tpi (white), and 
tp = (gray) emerged, but in the final stage t > 10^ the 
variant with tp = —tpi disappeared here. 



IV. NUMERICAL RESULTS 

We numerically study the dynamics of our model. We 
may demonstrate the validity of our equilibrium theory 
in steady states attained at long times. In our system tp is 
a nonconserved variable obeying the relaxation equation. 



d SF 



(4.1) 



where 5F/Sip is given in Eq.(2.14) and Lg is a con- 
stant. We integrated the above equation in 2D under 
the periodic boundary condition. We assume u > and 
(3/u = 1.5. Then, for h — 0, our theory predicts tpi — 
0.612Mo, Vo = 0, T,,/ro = 0.047, (r,,-T„)/ro = -0.516, 
Ko£ = 0.354, and ki£ = 0.596. @ These values will be 
compared with those from our simulations. 



A. Isotropic elasticity 

We first assume the isotropic elasticity. We measure t, 
h, and tp in units of tq, ho, and Mq, respectively, where 

To = uyv, ho = v{ulvf'^, Mo = {u/v)^/^ (4.2) 



W / M 



t=(iOOO 




FIG. 9: Time evolution of tp after changing r from —to to 
— 0.15ro at t = (upper panel) and steady profile of tp/Mo 
obtained at t = Qx 10^ (lower panel) for P/ti = 1.5 and h = 
in isotropic elasticity. In the phase ordering, ordered regions 
with tp = tpi (black) and disordered regions (gray) emerged. 
A circular ordered domain remained at long times in this run. 



Here M/Mq = [3(/3/?2 - 
space and time are 



l)/4]i/2 from Eq.(3.9). Units of 



to = LoTo, e = (C/to)1/2. 



(4.3) 



The scaled time to^t and the scaled space position £^^r 
are simply written as t and r to avoid cumbersome nota- 
tion. The system size is 200 x 200 and the mesh length 
is £, so the system length is 200^. In terms of the scaled 
order parameter ^' = tp/Mo, Eq.(4.1) is rewritten as 



'dt 



To U 



-A. (4.4) 

ho 



As the initial condition at i = 0, ^' at each lattice point 
consists of a homogeneous constant and a random num- 
ber in the range [—0.01,0.01]. 

In Fig. 8, we show the phase ordering process from 
a disordered state to a coexisting state. At t = 0, ^' 
was a random number. For t > we lowered t from 
to — 0.3ro to induce phase ordering. From our theory, 
this final t is in the coexisting window [tcx — t^, t^x] and 
the predicted average order parameter is 0.378Mo with 
(p = 0.617. Since h = and (tp) = at i = 0, the two 
variants with tp = ±tpi appeared in the early stage, but 
the ordered domains with tp = —tpi disappeared in this 
run when the domain size became of the order of the 
system size. (In other runs the variant with tp = tpi dis- 
appeared as well.) In the steady two phase coexistence 
at i = 10^ (lower panel in Fig. 8) interfaces are horizon- 
tal (parallel to the x axis), where tp = 0.612Afo in the 
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i=10000 




200 



FIG. 10: Steady profile of ip/Ma in two phase coexistence 
obtained aX t = 10* in isotropic elasticity, where (3/u — 1.5, 
r = — 0.13to, and h — O.Qhc- Here the system is close to the 
field-induced critical point and the interface region is broad- 
ened. 




t=250 



1 



t-JOOO 




ordered phase and {tp) = 0.397Mo. The former coincides 
with the predicted value, while the latter is slightly larger 
than predicted. 

In Fig. 9, we show the phase ordering process from 
a one phase state at r = — tq to a coexisting state at 
T = — 0.15ro at /i = 0. That is, at f = 0, 4* was the 
sum of the equilibrium one phase value 0.786 determined 
by Eq.(2.24) and a random number. The final t here 
is higher than the lower instability value — 0.313ro in 
Eq.(2.23). Hence phase ordering should take place into a 
coexisting state where (f> — 0.350Mo and (ip) — 0.214Mo 
are predicted. In the simulation, regions of the disordered 
phase appeared, while tp in the ordered phase changed 
to Tp = tpi. In the steady two phase coexistence at 
t = 6 X 10'^ (lower panel in Fig. 9), a circular ordered 
domain was realized. There, we find ijj — 0.594A/o in 
the domain and {tjj) — 0.232Aio. These values are only 
slightly different from those predicted. 

In Fig. 10, we present a steady profile oitp a.t h = 0.9/ic 
and T = — 0.13ro, where the system is close to the critical 
point in Eqs.(3.12)-(3.14) and the interface thickness is 
much widened. For P/u = 1.5 and at this field, our 
theory gives ^i/Mq = 0.461, Vo/Mq = 0.304, t^^/to = 
-0.032, {tc^ - t^)/tq = -0.212, Koi = 0.149, and kiI = 
0.164. For the t adopted, we predict cp — 0.543 and 
(ip) = 0.389A/o. In the simulation, the maximum and the 
minimum of ip are 0.457Afo and 0.302Mo, respectively. 
These values are very close to the theoretical values of 
tpi and V'o- Furthermore, the observed average {tp) = 
0.392Mo is also close to its theoretical average, though 
the interface regions are very wide here. 



B. Cubic elasticity 

Next we integrate Eq.(4.1) in 2D on a cell of 256 x 256 
assuming the cubic elasticity with Cn — C12 = C44 = K , 



FIG. 11: Time evolution of tp after changing r from to 
— 0.3ro at f = (upper panel) and final steady profile of ip/Mo 
obtained at t = 10* (lower panel) in cubic elasticity for pS/u — 
1.5 and ^ = 0. As in Fig. 8, three regions with tp = tpi (black), 
tp = —tpi (white), and tp = (gray) emerged in the initial 
stage. Interfaces tend to be parallel to the x ox y axis. 



where K = (Cn + Ci2)/2. Then = -1 from Eq.(2.11) 
and the softest directions are [10] and [01]. As in the 
isotropic case, space and time are measured in units of i 
and to in Eq.(4.3) and we set /3 = 2a^ /Cn = 1.5u > 0. 
The mesh size of integration is i. In terms of the scaled 
order parameter \1/ = tp/tpQ, the dynamic equation in the 
2D cubic case is written aal i^^i^^i^^ 



9* 
'dt 



Vlf2 _ vI/4 

To 



-. (4.5) 

ho 



From Eqs.(2.8) and (2.14) we express G(r) in the Fourier 
expansion. 



(4.6) 



where is the Fourier component of $ = 4"^ — (^'^) 

and w{k) = w{9) in Eq.(2.13) depends on the angle 
defined by cos6' = k^/k. 

In Fig. 11, we lowered r from to — 0.3ro at ft, = as 
in Fig. 8. Here the anisotropy of the domain structure 
arises from the angle dependence of w{0) in Eq.(4.6). 
In the steady state in the lower panel, the maximum of 
tp is O.6I3M0 and the average {tp) is 0.393A/o, in close 
agreement with the predicted values and those in Fig. 8. 

In Fig. 12, we show a steady profile of -0 for r = 
— 0.13ro and /i = as in Fig. 9. Here a square ordered 
domain is embedded in a disordered region in equilib- 
rium. In the figure, the maximum and the average of tp 
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t=10000 




250 



FIG. 12: Steady square profile ol tp/Mo in two phase coexis- 
tence obtained at t = 10* in cubic elasticity, where f3/u — 1.5, 
T — — 0.15to, and h = 0. 



b/h=os 




FIG. 13: Steady, one-dimensional curves of ip/Mo in two 
phase coexistence with f3/u = 1.5 near the field- induced criti- 
cal point in cubic elasticity. Here r — — O.lSro and — 0.15ro for 
the two curves of h — 0.9/ic, while r = —O.lSro for h — 0.95/ic. 
The more ordered region expands with lowering r at fixed h. 



are 0.586Mo and 0.260Mo, respectively. The former is 
slightly smaller than the predicted value 0.612Mo, while 
the latter is considerably larger than the predicted value 
0.214Mo. 

In Fig. 13, we show one-dimensional steady pro- 
files changing along the x axis near the critical point. 
The maximum, the minimum, and the average of ■0 are 
(0.461,0.305,0.392) for h/h^ = 0.9 and t/tq = -0.13, 
(0.461,0.308,0.409) for h/h^ = 0.9 and t/tq = -0.15, 
and (0.439,0.336,0.398) for h/hc = 0.95 and t/tq = 
—0.13. These values closely agree with those from our 
theory. In these one-dimensional cases, the profiles coin- 
cide with those in the isotropic case. 



V. SUMMARY AND CONCLUDING REMARKS 

We have examined the phase transition behavior of 
compressible Ising models at fixed volume in the mean 



field theory. In our model the order parameter ^ is 
isotropically coupled to the dilation strain ei as ip'^ei in 
the free energy, which is the simplest case. Nevertheless, 
complicated phase behavior follows at constant volume. 
We summarize our main results. 

(i) We have found two phase coexistence in a closed re- 
gion in the r-h plane as in Fig. 3. The coexistence region 
appears under the condition u < (3 given in Eq.(3.1). If 
u > and f3 is not large, it can be satisfied near the tri- 
critical point. If u < 0, it can occur even away from the 
tricritical point. 

(ii) The order parameter values in the two phases, ■01 
and ipOj a-rs determined by h only and is independent 
of T as in Fig. 2. The average order parameter {ip) — 

-I- (1 — 0)1/19 is increased smoothly as r is decreased 
in the window region Tcx — Tw < t < Tqx for h < 
since the volume fraction depends on t as in Eq.(3.16). 
The average order parameter (■0) and the susceptibility 
X — d{tp)/dh are displayed in Figs. 4-6. The specific heat 
Cv is a constant in two phase coexistence as in Eq.(3.26). 

(iii) At the field-induced critical point h = he and r = Tc, 
the correlation length 1/k grow and the surface tension 
7 goes to zero as in Fig. 7, while x does not diverge. 

(iv) We have integrated the dynamic equation, which is 
Eq.(4.4) for the isotropic elasticity and Eq.(4.5) for the 
cubic elasticity. A change of r from the one phase region 
into the unstable region induces phase ordering as illus- 
trated in Figs. 8-13. It can occur with decreasing t as 
in Figs. 8 and 11 and with increasing r as in Fig. 9. In 
the final two phase states, the values of ^0 and its space 
average closely agree with the theoretical values. 

We make some further remarks. 

(i) At constant pressure, two phase coexistence occurs 
only on a line in the r-h plane as in the rigid lattice 
case, but phase separation can be much affected by 
the elastic coupling (see the appendix)^. It is worth 
noting that the transition depends on the sample shape 
in hydrogen-metal systems at constant pressureii^^, 
where the proton concentration is linearly coupled to 
the dilation-^. 

(ii) We mention Monte Carlo simulations on a bi- 
nary alloy by Landau's group ^S^iiii^. They assumed 
that a mixture undergoing unmixing corresponds to 
fcrromagnets and that forming a superstructure to 
antifcrromagnets. In these cases, different results fol- 
lowed in the fixed volume and fixed pressure conditions. 
However, the unmixing transition in the presence of the 
size differencei is not isomorphic to the ferromagnetic 
transition. In the former the linear coupling^"^ appears 
between the concentration c and ei in the form "061, 
while in the latter the exchange interaction does not 
break the invariance of V' — * —'0 a-i^d the elastic coupling 
is quadratic as ip'^ei. At present we cannot compare our 
theory and their simulations. 

(iii) Yamada and Takakura numerically solved a 
time-dependent Ginzburg-Landau model for an order 
parameter and a strain in one dimension. They found 
appearance of a disordered region in a lamellar ordered 
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region-"^. Their finding is consistent with our theory. 

(iii) In real metamagnets, there is no field conjugate to 
the antiferromagnetic order and the tricriticality has 
been realized by changing magnetic field or hydrostatic 
pressure. At fixed volume, our theory predicts two 
phase coexistence in a temperature window near the 
tricritical point and near the line of first order phase 
transition. From Eq.(3.21) the width of the window 
sensitively depends on the coupling constant a as 
T^j Aq = 3l3{/3 — u)'^ /4:vAq, where Aq is the coefficient in 
Eq.(2.3) and /? = 2a^/K. 

(iv) In our mean field theory, we have neglected the 
renormalization effect near the critical point, which 
can be intriguing in the presence of the cubic elastic 
anisotropy®'^ . It should be further studied together with 
the influence of the global elastic constraint studied in 
this work. 

(v) We should generalize our theory to more complex 
systems. At the ferroelectric transitioni^,, the polar- 
ization vector is coupled to the strains. In binary 
alloys, phase separation and an order-disorder phase 
transition can take place simultaneously^, where the 
concentration c and the structural order parameter ^ 
are both coupled to ei in the form {aic + a2'4>'^)ei in the 
free energy^!. There can also be a number of anisotropic 
elastic couplings between the order parameter and the 
tetragonal or shear strain. We will soon report on phase 
transition including a Jahn- Teller coupling^^. 
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Appendix A: Fixed pressure condition 



We here eliminate the elastic field at fixed 
pressure^iiiii. Under isotropic appfied stress, we 
assume an isotropic average dilation change (ei) caused 
by the order parameter change. The average stress 
should be unchanged from that in the reference state, so 
we require (cry ) = in Eq.(2.6) to obtain 



-a(^')/A^ 



(A.l) 



in terms of the bulk modulus K . We impose the periodic 
boundary condition on the deviation, 5ui — Ui — {ei)xi/d, 
whose Fourier component can be expressed in terms of 
fj^ in the same form as that of Ui in the fixed volume 
case. The free energy consists of F in Eq.(2.15) and 

AF = -Va^{^l/f/2K. (A.2) 
The total free energy F' ^ F + AF is written as 



F' = dr 



(A.3) 



where f = ,fo — a^tp'^/2K and _B is a positive coefficient, 
B ^2a^/K -2wM- (A.4) 

Here wm is given by Eq.(2.17). The positivity of B arises 
from Cii — Ci2 > and C44 > 0. The one phase or- 
dered states are determined by /. The same form of the 
free energy was derived by Littlewood and Chandra^l for 
BaTiOa, who argued that the term proportional to B can 
much decrease the nucleation rate from the paraelectric 
to ferroelectric state. In our problem, we draw the fol- 
lowing conclusion in the mean field theory. In the fixed 
pressure condition, there can be two phase coexistence 
only on a first-order coexistence line in the r-h plane. In 
fact, (/) in Eq.(3.3) would be minimized for = or 1 
outside the coexistence curve if positive /3 were replaced 
by negative —B. 
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